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Abstract Complex network theory provides a powerful framework to statistically 
investigate the topology of local and non-local statistical interrelationships, i.e. 
teleconnections, in the climate system. Climate networks constructed from the 
same global climatological data set using the linear Pearson correlation coefficient 
or the nonlinear mutual information as a measure of dynamical similarity between 
regions, are compared systematically on local, mesoscopic and global topological 
scales. A high degree of similarity is observed on the local and mesoscopic topo- 
logical scales for surface air temperature fields taken from AOGCM and reanalysis 
data sets. We find larger differences on the global scale, particularly in the be- 
tweenness centrality field. The global scale view on climate networks obtained 
using mutual information offers promising new perspectives for detecting network 
structures based on nonlinear physical processes in the climate system. 

1 Introduction 

During the last decade, the development and application of complex network theory generated 
a wealth of novel insights into the nature of complex systems in various areas of science, e.g. 
the internet and world wide web in computer science, food webs, gene expression and neural 
networks in biology and citation networks in social science [1-3]. The intricate interplay be- 
tween the structure and dynamics of real world networks has received considerable attention [4] . 
Particularly, synchronization arising by the transfer of dynamical information in complex net- 
work topologies has been studied intensively [5] . The application of complex network theory to 
climate science is a very young field, where only few studies have been reported recently [6-13]. 
The vertices of a climate network are identified with the spatial grid points of an underlying 
global climate data set. Edges are added between pairs of vertices depending on the degree of 
statistical interdependence between the corresponding pairs of anomaly time series taken from 
the climate data set. 

When studying networks in the climate system, one has to assume that its dynamics can be 
approximated reasonably well by a grid of low dimensional nonlinear dynamical systems inter- 
acting only with their spatial neighbors according to the locality principle of classical physics. 
Note that this assumption is made implicitly, when the fundamental partial differential equa- 
tions of fluid mechanics are discretized and integrated in large scale climate simulations by the 
coupled atmosphere-ocean general circulation models (AOGCMs) used in climate science. Due 
to the continuity of the underlying physical fields, such as temperature or pressure, neighbor- 
ing grid points are dynamically correlated; these trivial local correlations usually decay quickly 
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within a typical length scale. Additionally, richly structured long range correlations appear, that 
were named teleconnections by the climatological community and have been studied extensively 
since the end of the nineteenth century [14]. 

The climate network approach enables novel insights into the topology and dynamics of the 
climate system over many spatial scales disclosed by local network measures, e.g. the number 
of first neighbors of a vertex v (the degree centrality k v ), mesoscopic measures such as the 
clustering coefficient and global measures, e.g. the average path length. The local degree cen- 
trality and related measures have been used to identify super-nodes (regions of high degree 
centrality) and to associate them with teleconnection patterns in the atmosphere, most notably 
the North Atlantic Oscillation (NAO) [6-8]. On the global scale, climate networks were found 
to possess 'small-world' properties due to long range connections (edges linking geographically 
very distant vertices), that stabilize the climate system and enhance the information transfer 
within it [6-8]. We stress, that the transfer of information in any complex physical system, e.g. 
the climate system studied here, will be carried by a flow of matter and energy. By studying 
the prevalence of long range connections in El Nino and La Nina climate networks [9] and the 
time dependence of the number of stable edges [10,11], it has been shown very recently, that 
the El Nino-Southern Oscillation (ENSO) has a strong impact on the stability of the climate 
system. 

In all works mentioned above, researchers have used the linear cross-correlation function of 
pairs of anomaly time scries to quantify the degree of statistical interdependence between differ- 
ent spatial regions. But the highly nonlinear processes at work in the climate system call for the 
application of nonlinear methods to obtain more reliable results. In a recent work on structures 
in the bctwecnness centrality field of climate networks [13], we have introduced mutual infor- 
mation [15] as a measure of statistical interdependence to climate network construction. The 
mutual information allows to capture nonlinear relationships between time series. We found 
that, while many properties of climate networks generated using the Pearson correlation and 
the mutual information at zero lag are qualitatively and quantitatively similar, the betweenness 
centrality field shows much greater deviations between the two construction methods. To check 
the possibility, that these pronounced differences are a signature of nonlinear processes in the 
climate system, and to bridge the gap between our nonlinear network construction method 
and the techniques previously used, we present a systematic statistical similarity study of the 
resulting climate networks. We show, that over a wide range of relevant edge densities (the 
fraction of the maximum number of possible edges present in the network), a high degree of 
similarity is maintained on local and mesoscopic topological scales. Furthermore, we address 
some of the more pronounced differences on the global topological scale, that are uncovered by 
bctwecnness centrality, and their possible relation to nonlinear processes in the climate system. 

The organization of the paper is the following: We first describe the data and the filtering 
and normalization procedures applied to it (Sect. [2]). After introducing the required elements of 
complex network theory (Sect.j3J), we proceed to develop in detail the method of climate network 
construction (Sect.|4|. In Sect.[5| we present the systematic comparison of the measures obtained 
from Pearson correlation and mutual information climate networks, respectively. Furthermore 
we provide a climatological interpretation. Some conclusions are drawn in Sect. [6] 



2 Data 

2.1 Description 

We utilize the monthly averaged global surface air temperature (SAT) field for climate network 
construction to maintain consistency with earlier works that analyzed the same field [6-11, 
13]. The SAT field allows to directly capture the complex dynamics on the interface between 
ocean and atmosphere due to heat exchange and other local processes. SAT therefore enables 
us to study atmospheric as well as oceanic dynamics within a common framework. We use 
reanalysis data provided by the National Center for Environmental Prediction/National Center 
for Atmospheric Research (NCEP/NCAR) [16] and model output from the World Climate 
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Table 1. Properties of global surface air temperature data sets. 



NCEP/NCAR rcanalysis | HadCM3 



Temporal coverage Jan 1948 - Dec 2007 Jan 1860 - Dec 1999 

T [months] 720 1680 

AX [°] 2.5 2.5 

A<t> [°] 2.5 3.75 

N 10224 6816 



Research Programme's (WCRP's) Coupled Model Intercomparison Project phase 3 (CMIP3) 
multi- model data set [17]. For optimal comparability with the reanalysis data, we choose a 20th 
century reference run (20c3m, as defined in the IPCC AR4) by the Hadley Centre HadCM3 
model. A data set consists of a regular spatiotemporal grid with time series Xi(t) associated to 
every spatial grid point i at latitude Aj and longitude fa. Start and end dates, length of time 
series T, latitudinal resolution AX, longitudinal resolution A(f> and the number of vertices of 
the corresponding global climate network N are given in Table [T] Note that we remove the 
polar grid points at A 6 {—90°, 90°} from the data sets, since the poles are represented by rows 
of grid points with identical dynamics. 

2.2 Filtering and normalization 

To minimize the bias introduced by the external solar forcing common to all time series in the 
data set, we calculate anomaly values, i.e. remove the mean annual cycle by phase averaging. 
Relabeling the time series by month to € {1, . . . , 12} and year y mapping Xi(t) — > ccj(y, to) one 
obtains anomaly time series at(y, to) = Xi(y, m) — (xi(y, m)) , that are consequently subjected to 
the inverse mapping a,(y, to) — > a,(t). Here and in the following (f(x)) x denotes the expectation 
value of observable / taken with respect to the variable x. Note that the anomaly time series 
already have zero mean. We furthermore normalize the anomaly time series to unit variance. 
Up to this point, we follow the method used previously by [9,10]. It is known, that the annual 
cycle induces higher order effects such as seasonal variability of anomaly time series variance. 
We find that using only data from a particular season to avoid biases due to this effect does not 
alter our results substantially, so that we choose to use the whole data set for a more accurate 
evaluation of statistical interdependence. 



3 Elements of complex network theory 

Formally, a network or graph is defined as an ordered pair G := (V, E) containing a set V = 
{1, N} of vertices together with a set E of edges {i, j}, which are 2-element subsets of V. In 
this work we consider undirected and unweighted simple graphs, where only one edge can exist 
between a pair of vertices and self-loops of the type {i,i} are not allowed. This type of graph 
can be represented by the symmetric adjacency matrix 



The edge density of a network is given by p = \E\/(f) = (k v ) v /N, \E\ being the number of 
edges in the graph and (k v ) v the mean vertex degree. The network measures defined below 
were selected for this study, because they allow us to compare different aspects of climate 
network topology on local, mesoscopic and global scales (Table [2| and are well established 
in the literature [2-4, 18]. Degree centrality, the related area weighted connectivity and the 
Hamming distance use only local information on the direct neighborhood of a vertex v. In 
contrast, the closeness and betweenness centralities as well as the average path length include 
global topological information by relying on shortest paths between pairs of vertices in the 
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Table 2. Classification of network measures into topological scales, fields and scalar measures. 





local 


mesoscopic 


global 


field 


k v , AWC V 


C v 




scalar 


H(A,B) 


C 


C 



network. This is why we refer to the latter three as global measures. On the mesoscopic scale, 
the local and average clustering coefficient depend only on information about neighbors and 
next neighbors of vertices. The concept of topological scales is elaborated in greater detail 
in [19]. We refer to measures assigning a real number j„ g R to each vertex v € V via a 
mapping V — > M : v i— » g v as fields. Scalar measures produce a single real number for the whole 
graph. 

Note that for the data sets analyzed here (Sect. [2]), vertices are not distributed homoge- 
neously on the earth's surface. The density of vertices increases from the equator towards the 
poles. This induces an inherent bias in the network measures studied, which prompts to use 
area weighted generalizations of the standard complex network measures, e.g. area weighted 
connectivity is the generalization of degree centrality. We have performed extensive studies of 
climate networks constructed from data interpolated to different grids and resolutions and find, 
that our results (Sect.pjf are not altered significantly by the vertex density bias [20]. This holds 
particularly for the highly interesting path based measures on the global topological scale. 



3.1 Local measures 



3.1.1 Degree centrality 



The degree or degree centrality [18] k v gives the number of first neighbors of a vertex v and 
can be calculated from the network adjacency matrix Ay using 

JV 

k v = S A v i. (2) 

i=i 

Vertices with exceptionally high degree centrality are usually referred to as hubs or super- 
nodes. We extend the use of this term to regions of spatially adjacent vertices with high degree 
centrality. 



3.1.2 Area weighted connectivity 

The area weighted connectivity 

AWC V = % 4 ^^ i (3) 

is closely related to the degree centrality k v of v. It corrects for the fact that in geographical 
networks defined on a grid, vertices correspond to regions of different area on the earth's surface. 
For the angularly equidistant grids considered in this work, the corresponding area of vertex 



v is proportional to the cosine of latitude A„ (see Sect. 2.11. AWC V can be interpreted as the 
fraction of the earth's surface area a vertex is connected to [7]. AWC is thus normalized to 
< AWC„ < 1. 
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3.1.3 Hamming distance 



The Hamming distance H(A, B) of two labeled simple graphs with adjacency matrices Ai 3 - and 
Bij measures the fraction of edges that have to be changed to transform one graph into the 
other [21]. Both graphs must contain the same number of vertices N. Specifically, H(A,B) is 
given by 

H(A,B) = {XOR(A ij ,B ij )) ij , (4) 

where 

XO« Wi , Ba) = {j ^* B " (5) 

Hamming distance is bounded by < H(A, B) < 1 and measures the global probability of 
non-equal entries in the two adjacency matrices. In our application we calculate the Hamming 
distance of two graphs with approximately equal edge density p. 

To evaluate the significance of this measurement, we compare it with the expected Hamming 
distance H R (p) of two independent Erdos-Renyi random graphs of edge density p [22]. The 
probability that the entries A i3 and B i3 differ between the two random graph adjacency matrices 
is given by p (Ay ^ B i3 ) = p {A i3 = 1) p (fly = 0) + p (A i3 = 0) p (fly = 1) = p(l - p) + (1 - 
p)p = 2p(l — p). Since all entries are independent, taking the expectation value reveals the 
expression H R (p) = (p (A i3 ^ fly))y = 2p(l — p). The expected Hamming distance H R (p) 
gives a reference point against which to judge the similarity of two graphs. We will make use 
of it in Sect. |5.3| to compare the performance of two network measures in climate network 
construction. 



3.2 Mesoscopic measures 



3.2.1 Local clustering coefficient 



We refer to C v as the local topological clustering coefficient or Watts-Strogatz clustering coef- 
ficient [1] of a vertex v. It gives the probability, that two randomly chosen first neighbors of v 
are also neighbors. With r v being the set of first neighbors of v and e(r v ) the number of edges 
connecting the vertices within the neighborhood F v , the clustering coefficient can be written as 



C (6) 

where the binomial coefficient ( ^) = ^k v (k v — 1) gives the maximum number of edges in r v . 
The local clustering coefficient is normalized to < C v < 1. 



3.2.2 Global clustering coefficient 



We speak of the (global) clustering coefficient C as the mean Watts-Strogatz clustering coeffi- 
cient 



C = (C v ) 



(7) 
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3.3 Global measures 

3.3.1 Closeness centrality 

Closeness centrality CC V measures the inverse average topological distance of vertex v to all 
others in the network [18], 

N - 1 

l^ii—l a vi 

where the topological distance or shortest path length dij is the minimum number of edges 
that have to be crossed to travel from vertex i to vertex j (d vv = by definition). If i and j 
arc not connected, the maximum topological distance in the graph d^ = N — 1 is used in the 
sum. Closeness centrality is normalized to < CC V < 1. Following our definition, CC V is large, 
when v is topologically close to the rest of the network. One should bear this in mind, because 
some researchers have used the inverse of our definition [18,19]. 

3.3.2 Betweenness centrality 

Assume that information travels through the network on shortest paths. There are <7jj shortest 
paths connecting two vertices i and j. We then regard a vertex v to be an important mediator 
for the information transport in the network, if it is traversed by a large number of all existing 
shortest paths. Mathematically, the betweenness BC V can be expressed by 

BC V =±^ (9) 

where <Jij(v) gives the number of shortest paths from i to j, that include v [18]. Here the 
contribution of shortest paths is weighted by their respective multiplicity atj . 

3.3.3 Average path length 

The average or characteristic path length £ of a graph is defined as the average topological 
distance between all pairs of vertices, 

£ = AE^- ( 10 ) 

\2) i<j 

Disconnected pairs of vertices are not included in the average, for a detailed discussion see [2]. 

4 Constructing climate networks 

To clarify the physical rational behind our method of climate network construction, we discuss 
it within the framework of synchronization from dynamical systems theory [23] . In a discretized 
model of the climate system, dynamical correlations can be envisioned as arising by (partial) 
synchronization of nonlinear oscillators on the grid that physically form a locally connected 
network. Even this simple network topology can generate nontrivial spatial patterns of synchro- 
nization [5,24,25]. The same is true for the synchronization of modes of variability in spatially 
continuous systems as the underlying fields of fluid- and thermodynamics [26], e.g. SAT. Many 
measures of synchronization have been proposed and used to infer coupling strength and di- 
rection between connected nonlinear oscillators [23,27]. The Pearson correlation coefficient [28] 
and the mutual information [29] were successfully employed to retrieve the network topology 
from the dynamics on the vertices alone. 
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The concept of synchronization provides a powerful paradigm to guide the enhancement of 
our understanding of the formation of (nonlinear) teleconnections in the climate system, and 
to stimulate the development of more advanced measures to detect these effects in measured 
data [23,26]. We hence propose that research aiming to construct networks from multivariate 
climatological data should be embedded within the framework of synchronization in complex 
networks [5]. 

4.1 Correlation measures 

In the spirit of simplicity facing comparably short time series and desiring consistency with 
the literature, we choose to first use the standard Pearson correlation coefficient and then 
cross-check the results by introducing mutual information to climate network construction. 
The mutual information will allow to investigate nonlinear dynamical relationships (nonlinear 
teleconnections) that are not fully detectable by using the linear Pearson correlation coefficient 
[30]. Note that we evaluate both measures at zero lag between time series. In principle, one 
can calculate a time delayed Pearson correlation (the cross correlation function) and mutual 
information [15]. This is appropriate when studying climate networks on smaller time scales 
using data sets with (sub-)diurnal resolution [10-12]. However, in the present work, we intend 
to study long term structural properties of the climate system on a scale of O(10 2 ) years using 
monthly averaged data. Most physical mechanisms of global information transfer in the SAT 
field, such as traveling Rossby waves, heat exchange between ocean and atmosphere or the 
advection of heat by surface currents in the ocean, act on time scales of less than one month. 
Therefore, it is reasonable to calculate the correlation measures at zero lag between anomaly 
time series. 



4.1.1 Pearson correlation coefficient 

The parametric empirical Pearson correlation coefficient Rij = (di(t)dj(t)) t — Rji estimates 
the strength of a linear relationship between two normalized time series &j and dj, given those 
are normally distributed. It produces spurious results for not normally distributed observables 
and nonlinear relationships. Consequently it should be used with care when constructing cli- 
mate networks. The non-parametric Spearman rank order correlation coefficient, that does not 
depend on the assumption of normally distributed observables, and are found to converge 
to the same value for nearly all pairs of time series taken from the data sets introduced in Sect. 
[2] The corresponding climate networks hence display close to identical network measures at all 
topological scales and we conclude, that utilizing the Pearson correlation coefficient to study 
linear climate networks is statistically justified here. 

In contrast to the standard definition of teleconnectivity [14] , we do not limit our analysis 
to strongly negative correlations. As in earlier works on climate networks, we use the absolute 
value of Pearson correlation P^j = \Rij\ = Pji to construct climate networks, since both large 
negative and positive values of Pearson correlation are indicative of a strong linear statistical 
interdependence . 

4.1.2 Mutual information 

In climate science, nonlinear measures of statistical interdependence have been successfully ap- 
plied to uncover strongly nonlinear relationships of climate observables, e.g. the phase coherence 
between ENSO and the Indian Monsoon [31]. Mutual information from information theory is 
another nonlinear measure now widely applied in many fields of science, ranging from linguis- 
tics [32] to computational neuroscience [29]. The mutual information Mjj can be interpreted 
as the excess amount of information generated by falsely assuming the two time series di and 
hj to be independent, and is able to detect nonlinear relationships [15]. By definition, My is 
large if the two time scries arc highly linearly (anti)correlated. In contrast, a strongly nonlinear 
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relationship between hi and hj yields large M^, but small Pij (see the upper left quadrant in 



Fig. 1(c) ). The mutual information can be estimated using 



where Pi((i) is the probability density function (PDF) of the time series hi, and Pij((J-, v) is the 
joint PDF of a pair (hi, hj). By definition, Mjj is symmetric, so that Mij = Mji. The standard 
unit of measurement of mutual information is the bit, if logarithms to base 2 are used. 

We use a simple histogram approach with equally siz ed bins for all pairs {i,j} to estimate 
the probability densities. Because the estimator (Eq. |11[ ) is known to depend on bin size and 
partitioning [33-35], we use an identical partitioning for all {i,j} to guarantee an optimal 
comparability of the Mij. We select a bin number of 64, that meets the Cochran criterion of at 
least 5 samples per bin for a typical time series length of 0(1 3 ). The basic algorithm applied 
here is computationally much less expensive than more advanced methods proposed in the 
literature [35,36], which is an important advantage when dealing with up to O(10 8 ) pairs in a 
global climate network. Our algorithm is feasible, since the application to network construction 
requires only the correct estimation of relative differences of Mij between all pairs of time series. 
In other words, in our application systematic under- or overestimation of mutual information 
is not a problem, as long as the error stays approximately constant across all pairs. 



4.2 Obtaining the network adjacency matrix 

We now construct the climate network by thresholding the correlation measure matrix Cij 
(Cij = Pij or C'ij = Mij), i.e. only pairs of vertices {i,j} that satisfy Cij > r are regarded as 
linked. By definition > 0, V{i,j} (see Sect, 
adjacency matrix Aij of the climate network is tf 



4.1|. Using the Heaviside function 0(x), the 



icn given by 

A ij = G{Ci j -r)-5 ij , (12) 

where <5y is the Kronecker delta. Note that inherits its symmetry from Cij and the resulting 
climate network is an undirected and unweighted simple graph. 



4.3 Choosing the threshold 

The last but nontrivial step in climate network construction is the selection of a threshold r, 
above which we consider a pair of vertices to be connected. From a statistical point of view 
it is desirable to only maintain connections that are statistically significant with respect to 
some reasonable test and reject those not meeting this criterion. Classical significance tests 
and randomization experiments have been used to assess the value of t for climate networks 
constructed using the Pearson correlation coefficient [6, 7, 9] . We build on these results testing 
against randomly shuffled time series, Fourier surrogates and twin surrogates [37]. Twin sur- 
rogates correspond to the null hypothesis of trajectories with random initial conditions on the 
attractor of the original time series and are found to give the strictest bounds on the significance 
of network connections detected using Pearson correlation and mutual information. 



4.3.1 On the role of teleconnections 

From the perspective of complex network theory, we intend to uncover interesting structures in 
the topology of the climate network. Different features of the underlying correlation measure 
matrix CV, will be revealed at different thresholds r. Consequently, the choice of r has to 
reflect a trade-off between the statistical significance of connections and the richness of network 
structures unveiled. For example, note the potentially interesting long distance edges with 
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Figure 1. (a,b) Frequency plot in the space of correlation measure dj and edge distance Uj = 
Rearth arccos(sin( Xi ) sin(Aj) + cos(Ai) cos(Aj) cos(0i — 4>i)) f° r au N(N — l)/2 = 23, 228, 928 pairs of 
time series in the global HadCM3 SAT data set. The apparent oscillations with edge distance are 
an artifact of the finite spatial resolution of the underlying grid, (c) Frequency plot in the space of 
Pearson correlation Py and mutual information Mtj. All plots are based on 2D-histograms with 10 4 
equally sized rectangular bins. The color bars indicate the common logarithm of frequency. Vertical and 
horizontal lines mark the thresholds corresponding to edge density p = 0.005 for P t j and My (Fig.j2|. 
The asterisk in (c) delineates the quadrant containing edges that exist in the mutual information, but 
not in the Pearson correlation network of p = 0.005, and hence are candidates for strongly nonlinear 
connections. 



high Pearson correlation an d mu tua l info rmation at edge distance I > 15000km in the global 
HadCM3 SAT data set (Fig. 1(a) and 1(b) I. They will only be included in the climate network, if 



the threshold r < 0.65 for the Pearson correlation network, or r < 0.3 in the case of the mutual 
information network. Long distance edges with high correlation measure or tclcconncctions 
are responsible for all interesting and non-trivial features of climate networks, such as small- 
world behavior, super-nodes or betweenness structures. Without them serving as spatial short 
cuts in the network, only the locally connected underlying grid remains. Ergo the inclusion of 
teleconnections must be a necessary criterion in the choice of the threshold in order to obtain 
interesting results in climate network analysis. 
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(a) 

Figure 2. PDF p(C) of the correlation measure matrices dj for the HadCM3 SAT data set. The 
vertical line indicates the threshold r yielding an edge density p(r) = 0.005, that is equal to the shaded 
area, (a) Pearson correlation, r = 0.682, b) mutual information, r = 0.398. 



4.3.2 Dependence of network measures on edge density 

Systematic studies show a smooth dependence of most climate network measures on r in the 
range of edge densities considered in this work. This implies that small uncertainties in the 
choice of the threshold will not lead to strongly deviating results within the complex network 
framework. Here we discuss the threshold dependence of edge density p(r), and the edge density 
dependence of clustering coefficient C(p), average path length C(p), number of components 
n c (p), relative giant component size S(p) and average relative non-giant component size (s(p)) 
(Fig. [3]). Here a component constitutes a maximally connected subset of vertices of the network, 
i.e. a connected subset of vertices that is not reachable from any other vertex in the network. 
The term giant component is usually reserved for the largest component containing nearly all 
of the vertices in the network [2]. S(p) in turn always measures the relative size of the largest 
component, even if its size becomes comparable to that of other components. 

The edge density p(r) decays approximately exponentially due to the shape of the PDF of 
the absolute value of the correlation measure p(C) (in the following we abbreviate Cij by C), 



P(r) 



dCp(C). 



(13) 



Note that p(r) is a monotonic decreasing function of r. Correlation measure distributions found 
empirically from climate data generally have a connected support (Fig. [2]), so that p(r) is strictly 
monotonic dec reasin g and induces a one to one correspondence between threshold r and edge 
density p (Fig. 3(a)). 



The clustering coefficient C is found to stay approximately constant at intermediate values 
of p and decays to ze ro for small p (Fig. 3(b)), w hen t he netwo rk decomposes into a larger 
number n c (Fig. |3(d"J| of smaller components (Fig. 3(e) and |3(f)| . The average path length L 
decays approximately as a power law with growing p and has discontinuities at edge densities 
p M , where r M = t(/? m ) equals the correlation measure CV, of edges {i,j} with a high edge 
betweenness ccntrality [2], i.e. that lie on many shortest paths between pairs of vertices (Fig. 
3(c) |. When t > t m , these shortest paths become considerably longer and components might 
decouple from the network's giant component. This effect leads to a de crease of £ for small p 
since the network decomposes into smaller disconnected components (Fig. |3(f)| and path lengths 
are measured only within the components. The formation of a giant component encompassing 
nearly all v ertice s at p ss 0.0012, where the giant component size increases from S ~ 0.5 to 
S ~ 1 (Fig. |3(e)[ ), goes along with discontinuities of C and (s(p)). Note that all vertices have 
joined the giant component at p 0.020 for the HadCM3 SAT Pearson correlation network 
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Figure 3. Network measures as a function of threshold and edge density for global HadCM3 SAT 
networks constructed using Pearson correlation, (a) Threshold dependence of edge density p(r), (b) 
edge density dependence of clustering coefficient C(p) and (c) average path length C(p). (d) Edge 
density dependence of the number of components n c (p), (e) giant component size S(p) and (f) average 
non-giant component size (s(p)). The vertical lines indicate edge densities of p — 0.005 and p = 0.01 
and corresponding thresholds. 
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(Fig. |3(d)| and at p k 0.028 for the corresponding mutual information network (not shown 
here) . 

At all edge densities considered in Sect.[5]the giant component size is of 0(1). The influence 
of the non-giant components on measures such as average path length and closeness centrality 
is therefore negligible in the regime studied here, since larger deviations are only expected for 
p < p. This range of edge densities in turn is not relevant for the conclusions drawn from 
the comparison presented in Sect. [5] To study this regime of very small edge densities in 
detail, measures more robust to disconnected components such as the local efficiency (related to 
closeness centrality) and global efficiency (related to average path length) should be considered 
[2]. We chose the definitions given in Sect. [3] to maintain consistency with the existing literature 
on climate networks. 



4.3.3 Pragmatic choice of r 

We think that the problem of selecting exactly the right threshold is not as severe as might 
be thought. Climate network analysis deals with topological properties of correlation measure 
matrices and aims at gaining new insights heeding this paradigm. In the climate system, it is 
furthermore not immediately evident which physical entities should take the role of vertices 
and edges in a complex network. This constitutes the main conceptional difference between our 
method and attempts of recovering an unknown physically existent network structure from ver- 
tex dynamics as in the study of the brain [28,29,38-40], where one can argue that a more natural 
identification of neurons and axons with the vertices and edges of a neural network exists. It is 
known that in the classical local description of geophysical fluid dynamics of atmosphere and 
oceans, i.e. the Navier-Stokes equations combined with thermodynamic equations, the network 
of physical interaction has the structure of a regular grid [41]. In a discretized model, the dy- 
namics at each grid point is only coupled to the grid points in the immediate neighborhood. The 
complex topology observed in climate networks should therefore be treated as a manifestation 
of structure formation, that allows for uncertainties in the choice of parameters such as r. 

In the spirit of the ideas elaborated in the above paragraphs, we choose to fix the edge den- 
sity p when comparing the properties of climate networks generated using different correlation 
measures. This will result in different thresholds r, because the empirical correlation measure 
distribution p(C) clearly differs between linear Pearson correlation and nonlinear mutual infor- 
mation (Fig.pl). The selection of p is in each case guided by the principle of balancing between 
structural richness and statistical significance outlined above. 



5 Results 

After having introduced our methodology for climate network construction, we proceed to the 
main aim of this study: A comparison of networks generated using the linear Pearson correlation 
coefficient and the nonlinear mutual information on local, mesoscopic and global topological 
scales. The edge density p is varied between p m ;„ = and p m ax — 0-1 in equally sized steps. 



Recall, that small edge densities correspond to high thresholds (Sec. 4.3 1. For increasing edge 
density, edges with decreasing correlation measure are added to the network. Consequently, 
climate networks with a very high edge density p > 0.1 are not expected to contain meaningful 
information for climate data analysis, because they contain many connections that are not 
statistically significant, i.e. that are much more likely to arise by chance. For example, Tsonis 
et al. use the Pearson correlation coefficient and a threshold of t = 0.5 in all of their works [6-9], 
which corresponds to an edge density of p 0.01 for the global HadCM3 SAT data set analyzed 
here. They report that according to the Student's t test, a value of Pij = 0.5 is statistically 
significant above the 99% level. In our recent work, we use an edge density of p = 0.005 [13]. 
This larger threshold corresponds to an even higher significance level, because it is less likely to 
be exceeded by the correlation measures calculated from pairs of one original and one surrogate 
time scries. 
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Table 3. Spearman's Rho r a (p) of area weighted connectivity (AWC), local clustering coeffcient (C), 
closeness centrality (CC) and betweenness centrality (BC) fields and Hamming distances H(p) and 
H R (p) calculated from Pearson correlation and mutual information networks at edge densities p = 0.005 
and p = 0.01 for the global HadCM3 SAT data set. 





p = 0.005 


p= 0.01 


r AWC 
1 s 


0.95 


0.88 


r c 

' s 


0.80 


0.81 


r cc 

' s 


0.98 


0.95 


r BC 
' s 


0.70 


0.59 


H 


0.001 


0.003 


H R 


0.010 


0.02 



We compare the properties of the complex networks obtained at each edge density level on 
local, mesoscopic and global topological scales. We enable a qualitative discussion of similarity 
by plotting the fields of area weighted connectivity (Fig. [4}, local clustering coefficient (Fig. 
[5]), closeness (Fig. [6]) and betweenness centrality (Fig. [7]) on a world map at fixed edge density 
p = 0.005. The local deviations of these fields calculated for Pearson correlation and mutual 
information climate networks are highlighted by normalized difference fields (Fig. [8]). For a 
quantitative comparison at all edge densities considered, we calculated the Spearman rank 
order correlation coefficient or Spearman's Rho r s (p) of the corresponding fields taken from the 
Pearson correlation and mutual information networks (Fig. 9(d)|and Fig 



10(d)). 



We chose to 

use the Spearman's Rho instead of the Pearson correlation coefficient for this task, because it is 
known to be more reliable when applied to data with non-Gaussian PDF. This is an important 
property, considering that some of the fields we are interested in have a highly non-normal 
frequency distribution (Sect 



5.1 and Sect. 5.3 1. Furthermore at each edge density step, we 



consider the Hamming distance between the networks on the local topological scale, whereas 
on the mesoscopic and global scale we compare global clustering coefficient and average path 
length. 

In the following we will illustrate the comparison for the HadCM3 SAT data set in detail 
5.3 and Fig. [4]j5l [6] 7j [8] [9]) . Only the quantitative comparison is presented for 



5.2 



(Sect.l5T| 

the NCEP/NCAR reanalysis SAT data set (Fig.|10|, since we are lead to the same conclusions 
as for the model data set. Finally we present climatological interpretations of the observed 



network structures (Sect. 5.4) 



5.1 Local comparison 

On the local topological scale, we find that Pearson correlation and mutual information climate 
networks are very similar at low edge densities. At p = 0.005, the area weighted connectivity 
(Fig. HI) field show s only small deviations by visual inspection, that are most pronounced in the 



tropics (Fig. 8(a)). The rank order correlation coefficient r* c reaches a maximum between 
p = 0.005 and p — 0.01 and decays for larger edge densities (Fig. |9(d"J| . We obtain high values 
for p — 0.005 and p — 0.01 (Table[3]). Note that for the climate networks studied, area weighted 
connectivity has a fat tailed PDF^7] . 

The Hamming distance H(p) is alway s sma ller than the expected distance H R (p) of two 
random networks at edge density p (Fig. 9(a)[ ). It is notable, that H(p) seems to go to zero 



tangentially to the p-axis, i.e. H'(p)\ p= o 0, whereas H R (p)\ p= o = 2. Therefore most of the 
edges with the highest Pearson correlation and mutual information values must coincide. From 
analytical considerations and Monte-Carlo simulations we find that the standard deviation of 
the PDF of Hamming distance between the two random networks is of 0(N^ 1 ) for N ^> 1. 
This means that the expected deviations from the mean H R (p) are of O(10" 4 ) for the climate 
networks considered here. The difference between measured Hamming distance and H R (p) is 
by one order of magnitude larger than these expected deviations (Table |3| . We hence conclude 
that the observed similarity of Pearson correlation and mutual information networks can be 
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considered statistically significant, with respect to the null hypothesis of random networks of 
the same size N, at all edge densities considered. Particularly, the results elaborated in this 
section show, that at the edge densities used in earlier works on climate networks [6-9, 13], 
Pearson correlation and mutual information give very similar results on the local topological 
scale. 



5.2 Mesoscopic comparison 

The local and global clustering coefficients also reveal a high degree of similarity on the meso- 
scopic topological scale. Analogously to AWC, the local clustering coefficient fields are nearly 
indist ingui shable (Fig. [5]). However, the largest deviations appear to cluster along coastlines 
(Fig. |8(b) |. This interesting finding can be understood by considering the qualitatively dif- 



ferent aynamics of SAT over oceans and continents, e.g. the on average much larger seasonal 
variability over continents. Along coastlines, the correlation length of the SAT field is thus 
smaller than that expected over continents or the ocean away from the coast. Hence Pearson 
correlation and mutual information have a higher probability to disagree on the existence of 
edges between spatially adjacent vertices (local edges) along the coastline. These local and 
mesoscopic deviations in network structure are detected by the local correlation coefficient C v , 



that is by design particularly sensitive on the mesoscopic topological scale (Sect. 5.4 1. 

The rank order correlation coefficient re aches a maximum between p — 0.005 and p = 0.01 
and decays for larger edge densities (Fig. |9(d)[ ). We obtain high values for p = 0.005 and 
p = 0.01 (Table |3]). The global clus tering coefficients show only small deviations of O(10~ 2 ) 
at all edge densities considered (Fig. |9(b)l We get C p (0.005) = 0.682, C M (0.005) = 0.678 and 



C p (0.01) = 0.657, C M (0.01) = 0.668rTne local clustering coefficient field is close to normally 
distributed. 



5.3 Global comparison 

We observe more interesting behavior at the global topological scale. Closeness centrality at 
p = 0.005 does not deviate much qualitatively and quantitatively across the two types of net- 
works considered (Fig. 6j, the largest differences are detected in the tropics with a tendency 



to decrease with latitude towards the poles, and most notably over South America (Fig. 8(c) I 
The bctwccnness centrality field shows more pronounced qualitative regional differences (Fig. 
[7]). For example, note the differing high betweenness structures over t he oceans, particularly 
over the East Pacific, the North Atlantic and arctic regions (Fig. [8(d)T ). The rank order corre- 
lation coefficients rf c and rf c decay more quickly than the ones on the local and mesoscopic 
topological scale and fluctuate around values of r^ c w 0.1 and rf c ~ 0.4 for larger edge den- 
ier than the Spearman's Rho 
1 1. Confirming earlier studies, 
he closeness field is normally 



sities (Fig. 9(d) |. At p = 0.005 and p — 0.01, rf c is notably smaller than the Spearman's Rho 
of the other fields considered, while rf c is close to unity (Table [3 1. Confirming earlier studies, 
we find that bctwccnness follows a fat tailed PDF [42] , whereas t : 
distributed. 

These results indicate, that betweenness centrality may quantify the local differences be- 
tween networks constructed using Pearson correlation and mutual information at the global 
topological scale, that could be traces of nonlinear physical processes in the climate system. 
That the greatest deviations are found between the betweenness centrality fields is plausible, 
because betweenness is by definition a very sensitive measure and can locally depend heavily 
on the existence or non-existence of a small number of edges in the network [43] . Consider for 
example a small set of edges, that are the only connections between two large communities in a 
network. The vertices on either end of these edges have a high betweenness centrality, because 
all shortest paths between the two communities must contain them. If the bridging edges are 
removed, the betweenness centrality of the beachhead vertices must decrease significantly, since 
they can now only participate in shortest paths within their own community. This sensitiv- 
ity of betweenness leads to a large dynamic range of 20 orders of magnitude for the global 
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Figure 4. Area weighted connectivity fields for global HadCM3 SAT networks at p — 0.005 (linear 
color scale) obtained using a) Pearson correlation, b) mutual information. The rank order correlation 
between the two fields is rf wc (0.005) = 0.95. 
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Figure 5. Local Watts-Strogatz clustering coefficient fields for global HadCM3 SAT networks at 
p = 0.005 (linear color scale) obtained using a) Pearson correlation, b) mutual information. The rank 
order correlation between the two fields is (0.005) = 0.81. 



Will be inserted by the editor 



17 




0° 45°E 90°E 135°E 180° 135°W 90°W 45°W 



.00912 .00923 .00935 .00946 .00958 .00969 .00981 .00993 .01004 .01016 .01027 .01039 .01051 .01062 

(a) 




0° 45°E 90°E 135°E 180° 135°W 90°W 45°W 

Closeness 




.01116 .01133 .01149 .01166 .01182 .01199 .01215 .01232 .01248 .01265 .01281 .01298 .01314 .01331 

(b) 

Figure 6. Closeness centrality field for global HadCM3 SAT networks at p — 0.005 (linear color scale) 
obtained using a) Pearson correlation, b) mutual information. The rank order correlation between 
the two fields is rf c (0.005) = 0.98. The white regions on the map correspond to vertices that are 
disconnected from the network's giant component. 
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0° 45°E 90°E 135°E 180° 135°W 90°W 45°W 

Betweenness (ln(BC + 1)) 



7 



1.5 3.01 4.51 6.01 7.52 9.02 10.52 12.03 13.53 15.03 16.54 18.04 19.54 21.05 

(a) 




0° 45°E 90°E 135°E 180° 135°W 90°W 45°W 

Betweenness (ln(BC + 1)) 



1.53 3.06 4.59 6.13 7.66 9.19 10.72 12.25 13.78 15.31 16.85 18.38 19.91 21.44 

(b) 

Figure 7. Betweenness centrality fields for global HadCM3 SAT networks at p — 0.005 (logarithmic 
color scale) obtained using a) Pearson correlation, b) mutual information. The rank order correlation 
between the two fields is rf c (0.005) = 0.70. 
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Figure 8. Normalized difference fields Ag v = — |/ \J {gZ) w {9m) w °f network measure fields g^ 
and g^, 1 , calculated from Pearson correlation and mutual information HadCM3 SAT climate networks, 
(a) Area weighted connectivity, (b) local clustering coefficient, (c) closeness and (d) betweenness. 



HadCM3 SAT network, that calls for a logarithmic scale to properly visualize the betweenness 
distribution (Fig. [7]) . 

The average path length (Fig 
C p (0.005) = 13.4, C M (0.005) 



9(c) I agrees closely, with deviations of O(10 1 ). We obtain 
13.5 and £ p (0.01) 



5, £ M (0.01) 



5.4 Climatological interpretation 



We give brief climatological interpretations of the network properties unveiled by our approach, 
since the main aim of this study is the comparison of linear and nonlinear climate network 
construction methods (Sect. [4]). Super-nodes found in the AWC field (Fig. [4} over the tropics 
and locally the mid-latitudes, were shown to be related to major atmospheric teleconnection 
patterns [8]. For example, the region of increased AWC in the North East Pacific is associated 
to the well-known Pacific North- American (PNA) pattern [14]. The El Nino cold tongue in the 
tropical East Pacific is clearly visible in the AWC field, as well as in all other fields considered 
(Fig.|5j[6]and0). 

The local clustering coefficient is found to be of 0(1) in a connected region in the equatorial 
Pacific as well as locally along continental coastlines (Fig. I5J. The former indicates a high 
degree of dynamical similarity in the tropical Pacific [7,8], that is possibly related to ENSO. 
The latter are more likely to be a signature of our climate network construction method along 
the coastline and visible on the mesoscopic scale only, that we discuss in Sect. |5.2| 

The contouring of the closeness field (Fig. [6]) nicely shows the latitudinally growing influ- 
ence of the Coriolis force. Pressure gradient forces are balanced by the Coriolis force in the 
mid-latitudes for large scale atmospheric flows. This balance vanishes in the tropics, because 
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(c) 



(d) 



Figure 9. Results for the quantitative comparison of Pearson correlation and mutual information 
climate networks for the global HadCM3 SAT data set, shown as a function of edge density p (100 
edge density steps), (a) Hamming distance H(p) (continuous line) and expected random Hamming 
distance H R (p ) (da shed line) between the two networks. The expected deviations from H R (p) are of 
O(10 -4 ) (Sect. 5.1 1, (b) Global clustering coefficient C F (p) of the Pearson correlation (continuous line) 
and C M (p) of the mutual information network (dashed line), (c) Average path length C F (p) of the 
Pearson correlation (continuous line) and C M (p) of the mutual information network (dashed line), (d) 
Spearman rank order correlation coefficients rf wc (p) for the area weighted connectivity (continuous 
line), 7"s (p) for the local clustering coefficient (crosses), r^ c \p) for the closeness centrality (dash-dotted 
line) and r BC (p) for the betweenness centrality fields (dotted line). 
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Figure 10. This figure shows the same statistics as Fig. [9] but evaluated for the global NCEP/NCAR 
reanalysis SAT data set. 



the Coriolis force decays as sin(A) when latitude A approaches the equator. The closeness field 
also shows that the tropics form the center of the SAT climate network, the associated ver- 
tices being topologically closer to the rest of the network than vertices in the mid-latitudes 
and arctic regions. This finding can be explained by considering the comparably regular dy- 
namics of the tropical SAT field leading to many edges between tropical vertices and the more 
irregular dynamics in the mid-latitudes and arctic regions that results in fewer edges within 
the mid-latitudes and arctic as well as between these regions and the tropics [44]. In a global 
climate network, it is hence more probable to find shorter shortest paths starting from tropi- 
cal vertices, while shortest paths originating in mid-latitude and arctic vertices are on average 
longer. Moreover, we point out the lower closeness over Australia and Greenland indicating 
that these land-masses also form pronounced clusters in the SAT climate network, even though 
the local clustering coefficient field shows that they are not as highly locally interconnected as 
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the equatorial Pacific. These differences in local connectedness among the detected dyna mica l 
clusters are caused by the qualitatively different dynamics over land and oceans (Sect. 5.2 1. 
The land-sea difference is globally detected by closeness centrality and AWC: Vertices over 
land masses are found to be on average less well connected and topologically more remote than 
those over the oceans. 

We observe highly localized linear structures in the betweenness field (Fig. [7]), some of 
which appear to resemble major surface ocean currents such as the California and Peru currents 
following the western coastline of the Americas, or the East Greenland, Norwegian and Canary 
currents. Note that some of these curr ent resembling structures are particularly visible in the 
betweenness difference field (Sect. 5.3 1, indicating that nonlinear processes might be involved 
in the formation of some of the structures. In analogy to the major communication channels of 
the internet, we refer to these betweenness structures as the backbone of the climate network, 
because a large fraction of the dynamical information exchanged via topologically shortest 
paths between all possible pairs of vertices {i, j} must pass the high betweenness regions. This 
is particularly true for information transported by advective processes, where the assumption 
of information traveling on shortest paths can be substantiated by extremalization principles. 
In our recent work we report the discovery of the backbone and its possible role in stabilizing 
the climate system [13]. 

Note that the region very close to the equator in the tropical East Pacific has a comparatively 
low AWC, closeness and betweenness, but a high local clustering coefficient. This indicates that 
this region forms a internally densely connected cluster in a network sense, i.e. it is dynamically 
highly interrelated but nearly detached from the rest of the network. We interpret it as a 
pronounced manifestation of the equatorial Coriolis barrier [41], that can also be observed 
weakly over the equatorial Indian and Atlantic Oceans. 

In agreement with [6-8] we find that Pearson correlation and mutual information climate 
networks possess properties of 'small-world' networks [1,45], i.e. a small average path length 
C -C N and a large clustering coefficient of 0(1) (Table [3j Fig. 9] and 10 1. Complex 'small- 
world' networks with comparable global properties are frequently found in nature, e.g. the 
internet, power grids, social and neural networks, and constitute the subject of study of an 
equally diverse collection of sciences. The small average path length can be explained by the 
influence of teleconnections. This indicates that perturbations of the regional dynamics (vertex 
dynamics) can on average quickly affect the whole globe via paths consisting of statistically 
highly interrelated pairs of regions (edges). It has been argued that this serves to stabilize the 
climate system and to enhance the information transfer within it [6-9]. If the climate network 
was only locally connected, in other words if all teleconnections were removed from it, the 
average path length would be of O(N) as that of a regular grid. The high clustering coefficient 
is due to the spatial continuity of the underlying physical fields (e.g. SAT), that leads to a 
prevalence of local triangles [46] . 



6 Conclusions 

In summary, we have performed a systematic study of the similarity of climate networks con- 
structed using the linear Pearson correlation and the nonlinear mutual information across local, 
mesoscopic and global topological scales. First, we have motivated the comparison of the two 
types of networks at equal edge densities. We have considered only low edge densities, that 
were shown to yield networks containing statistically highly significant edges as established on 
the basis of various significance tests. It has been then consistently shown for AOGCM and 
reanalysis surface air temperature data, that the networks agree well on the local and meso- 
scopic topological scales. Using the surface pressure field to construct climate networks also 
yielded qualitatively similar results and identical conclusions on these scales. For the surface 
air temperature data sets, we have found some interesting qualitative and quantitative devia- 
tions at the global scale using betweenness centrality. Even though there still is a high degree of 
similarity, the deviations arc highly localized and structured pointing at a possible involvement 
of nonlinear processes in their formation. 
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This work also demonstrates, that our method of calculating mutual information for rel- 
atively short time series is reliable at least for the strongly linear interrelations detected by 
the Pearson correlation coefficient. The global topological scale is of particular interest, since 
it opens novel perspectives for the understanding of climatological phenomena. For example, 
as applied to the climate networks discussed in this article, betweenness centrality allows to 
measure the importance of localized regions on the earth's surface for the transport of dy- 
namical information within a climatological field in the long term mean [13]. Further work is 
needed to establish, whether the observed deviations on the global topological scale could be 
due to nonlinear physical processes in the climate system, that are only detectable using mu- 
tual information. In the future, we plan to assess this problem by constructing climate networks 
using a novel method based on statistical significance, i.e. by adding edges to the climate net- 
work depending on the significance level of the correlation measure with respect to reasonable 
null hypotheses. One could then identify candidates for nonlinear interrelationships as edges 
that have an associated significant mutual information and a Pearson correlation that is not 
significant. 
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